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We present an unquenched calculation of the quark propagator in Landau gauge with 2+1 flavors 
of dynamical quarks. We use configurations generated with an improved staggered ( "Asqtad" ) action 
by the MILC collaboration. This quark action has been seen to have excellent rotational symmetry 
and scaling properties in the quenched quark propagator. Quenched and dynamical calculations 
are performed on a 20 x 64 lattice with a nominal lattice spacing of a — 0.125 fm. The matched 
quenched and dynamical lattices allow us to investigate the relatively subtle sea quark effects, and 
even in the quenched case the physical volume of these lattices gives access to lower momenta than 
our previous study. We calculate the quark mass function and renormalization function for a variety 
of valence and sea quark masses. 

PACS numbers: 12.38.Gc 11.15.Ha 12.38.Aw 14.65.-q 



I. INTRODUCTION 

Quantum Chromodynamics is widely accepted as the 
correct theory of the strong interactions and the quark 
propagator is its most basic quantity. In the low momen- 
tum region it exhibits dynamical chiral symmetry break- 
ing (which cannot be derived from perturbation theory) 
and at high momentum it can be used to extract the 
running quark mass Q. In lattice QCD, quark propa- 
gators are tied together to calculate hadron masses and 
other properties. Lattice gauge theory provides a way to 
calculate the quark propagator directly, providing access 
to quantities such as operator product expansion (OPE) 
condensates 2]. In turn, such a calculation can provide 
technical insight into lattice gauge theory simulations. 

The systematic study of the quark propagator on the 
lattice has also provided fruitful interaction with other 
approaches to hadron physics, such as instanton phe- 
nomenology 0, chiral quark models and Dyson- 
Schwinger equation studies 0, |(J. The lattice is a 
first principles approach and has provided valuable con- 
straints for model builders. In turn, such alternative 
methods can provide feedback on regions that are dif- 
ficult to access directly on the lattice, such as the deep 
infrared and chiral limits. 

The quark propagator has previously been studied us- 
ing Clover 0, @ , staggered Wi Holl and Overlap 0, 
actions. For a review, see Ref. [L|. All these actions have 
different systematic errors and the combination of these 
studies has given us an excellent handle on the possible 
lattice artifacts. In every case, however, they have been 
performed in the quenched approximation and have been 
restricted to modest physical volumes. 

In this paper we report first results for the quark prop- 
agator including dynamical quark effects. We use config- 
urations generated by the MILC Collaboration 0] avail- 
able from the Gauge Connection 0. These use "Asq- 



tad", 0(a 2 ) improved staggered quarks 16], giving us 
access to relatively light sea quarks. In the quenched ap- 
proximation, the quark propagator for this action has ex- 
cellent rotational symmetry and is well behaved at large 
momenta [jj. We use quenched and dynamical config- 
urations at the same lattice spacing and volume, which 
enables us to observe the relatively subtle effects of un- 
quenching. These lattices are also somewhat larger than 
those of previous studies, giving us access to smaller mo- 
menta. 



II. DETAILS OF THE CALCULATION 

The quark propagator is gauge dependent and we work 
in the Landau gauge for ease of comparison with other 
studies. Landau gauge is a smooth gauge that preserves 
the Lorentz invariance of the theory, so it is a popular 
choice. It will be interesting to repeat this calculation 
for the Gribov-copy free Laplacian gauge, and this is left 
for a future study. 

The MILC configurations were generated with 
the 0(a 2 ) one- loop Symanzik- improved Liischer-Weisz 
gauge action [Tt) . The dynamical configurations use the 
Asqtad quark action, an 0(a 2 ) Symanzik- improved stag- 
gered fermion action. They have two degenerate light 
fermions, for the u and d quarks, and a heavier one 
for the strange quark. We explore a variety of light 
quark masses, with the bare strange quark mass fixed 
at ma = 0.05, or m = 79 MeV for a = 0.125 fm [3. 
In all cases the Asqtad action is also used for the va- 
lence quarks. The values of the coupling and the bare 
sea-quark masses are matched such that the lattice spac- 
ing is held constant. This means that all systematics are 
fixed; the only variable is the addition of quark loops. 
The simulation parameters are summarized in Table [I] 

On the lattice, the bare propagator S(a;p 2 ) is related 



2 



to the renormalized propagator S Ten (/i;p 2 ) through the 
renormalization constant 



S(a;p 2 ) = Z 2 (a^)S ren (^p 2 ). 



(1) 



In the continuum limit, Lorentz invariance allows one to 
decompose the full quark propagator into Dirac vector 
and scalar pieces 

S- 1 (p 2 )=iA(p 2 )yp + B(p 2 ) (2) 

or, alternatively, 

S- 1 (p 2 )=Z- 1 (p 2 )[ il -p + M(p 2 )}, (3) 

where M(p 2 ) and Z{p 2 ) are the nonperturbative mass 
and wave function renormalization functions, respec- 
tively. Asymptotic freedom implies that, as p 2 — > oo, 
S(p 2 ) reduces to the free propagator 



s-Hp 2 



i"f ■ p + m, 



(4) 



up to logarithmic corrections. The mass function M is 
renormalization point independent and for Z we choose 
throughout this work the renormalization point as 3 GcV. 

The tree-level quark propagator with the Asqtad ac- 
tion has the form 



s-Hp) 



(5) 



where q{Pfj) is the kinematic momentum given by 

1 



9m 



sin(p M )[l + ^sin 2 ^) 



(6) 



The 7^ form a staggered Dirac algebra. Having identi- 
fied the kinematic momentum, we define the mass and 
renormalization functions by 

S~ l (p) = Z-\q) [i J2 (%)QM + M(q)] . (7) 



Complete details of the extraction of the mass and renor- 
malization functions from the Asqtad propagator are de- 
scribed in the appendix. 



TABLE I: Lattice parameters used in this study. The dynam- 
ical configurations each have two degenerate light (up/down) 
quarks and a heavier (strange) quark. The lattice spacing 
is a — 0.125(3) fm, where the uncertainty reflects the varia- 
tion of a over the set of lattices considered in this analysis. 
Bare light quark masses ma — 0.01, 0.02, 0.03, 0.04 correspond 
to masses of 16 — 63 MeV. The bare strange quark mass is 



ma 


= 0.05 


or 79 MeV. 








Dimensions 


P 


Bare Quark Mass 


# Config 


1 


20 3 


x 64 


8.00 


quenched 


265 


2 


20 3 


x 64 


6.76 


16 MeV, 79 MeV 


203 


3 


20 3 


x 64 


6.79 


32 MeV, 79 MeV 


249 


4 


20 3 


x 64 


6.81 


47 MeV, 79 MeV 


268 


5 


20 3 


x 64 


6.83 


63 MeV, 79 MeV 


318 



III. QUENCHED RESULTS 

First we compare our quenched results to some previ- 
ously published data obtained on a smaller lattice [lCj . 
All the data illustrated in the following are cylinder 
cut 0]. This removes points most susceptible to ro- 
tational symmetry breaking, making the data easier to 
interpret. As is well known, the definition of lattice spac- 
ing in a quenched calculation is somewhat arbitrary, and 
indeed the quoted estimate for our smaller ensemble is 
not consistent with that published for the MILC con- 
figurations. We determined a consistent value of the lat- 
tice spacing by matching the gluon propagator calculated 
on the old ensemble to that of the new ensemble [20| . 
This procedure yields a new nominal lattice spacing of 
a = 0.105 fm and physical volume of 1.7 3 x 3.4 fm for 
the old lattices. Examining the quark propagator on the 
two quenched ensembles, shown in Fig. ^ we see that 
the agreement is excellent. This indicates that both fi- 
nite volume and discretization effects are small. The flat- 
tening in the deep infrared of both scalar functions is a 
long-standing prediction of DSE studies 

We show results for the larger quenched lattice for a 
variety of bare quark masses in Fig. [21 Once again we see 
that for quark masses less than or approximately equal 
to that of the strange quark, the lowest momentum point 
of the mass function is insensitive to quark mass. 



IV. EFFECTS OF DYNAMICAL QUARKS 

Here we compare the scalar functions for the quenched 
and dynamical propagators. For a given bare mass, the 
running mass depends upon both the number of dynam- 
ical quark flavors and their masses. To make the most 
appropriate comparison we select a bare quark mass for 
the quenched case (ma — 0.01) and interpolate the dy- 
namical mass function so that it agrees with the quenched 
result at the renormalization point, q = 3 GeV. The re- 
sults are shown in Fig. [21 

The dynamical case does not differ greatly from the 
quenched case. For the renormalization functions, there 
is no discernible difference between the quenched and un- 
quenched cases. However the mass functions do reveal 
the effects of dynamical quarks. Dynamical mass gener- 
ation is suppressed in the presence of dynamical quarks 
relative to that observed in the quenched case. This is in 
accord with expectations as the dynamical quark loops 
act to screen the strong interaction. 

Further comparisons can be made in the chiral limit. 
In Fig. 01 both quenched and dynamical data have been 
extrapolated to zero bare quark mass by a fit linear in 
the quark mass. In the dynamical case, the extrapola- 
tion is limited to the case where the valence and sea quark 
masses match for the light quarks. While the results are 
qualitatively similar to those of Fig. [SJ increased separa- 
tion is observed in the generation of dynamical mass. As 
discussed above, for a given bare quark mass, the running 
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FIG. 1: Comparison of quenched wave-function renormaliza- 
tion and mass functions at approximately the same bare quark 
mass. The quark propagator from the 20 3 x 64 lattice with 
lattice spacing a = 0.125 fm at m — 47 MeV (open circles) 
is compared with the previously published quark propagator 
from a 16 3 x 32 lattice with lattice spacing a — 0.105 fm at 
m = 45 MeV (full triangles). The renormalization point for 
Z(q 2 ) is set at q = 3 GeV. 



FIG. 2: The quenched renormalization function (top) and 
mass function (bottom) for a selection of quark masses, in- 
cluding ma = 0.100, about twice the strange quark mass. 
The renormalization point for Z(q 2 ) is set at q — 3 GeV. 



mass is larger in full QCD than in quenched QCD. 

Fig. shows the mass and renormalization functions 
in the dynamical case for a variety of quark masses. Here 
the valence quark masses and the light sea quark masses 
are matched. The results show that the renormalization 
function is insensitive to the bare quark masses studied 
here. The results for the mass function are ordered as 
expected with the larger bare quark masses, m, giving 
rise to a larger mass function. 

Finally, we comment on the approach to the chiral 
limit. In Fig. we show the mass function for five dif- 
ferent momenta plotted as a function of the bare quark 
mass. The momenta considered include the lowest mo- 
mentum of 0.155 GeV and 0.310, 0.495, 0.700 and 0.993 
GeV to explore momentum dependent changes in the ap- 
proach to the chiral limit. At larger momenta, the mass 
function is observed to be proportional to the bare quark 
mass. However, at small momenta, nonperturbative ef- 



fects make this dependence more complicated. For exam- 
ple, a recent Dyson-Schwinger study predicts a downward 
turn as the bare mass approaches zero Q . 

For the lowest momentum points, nonlinear behavior 
is indeed observed. For the quenched case, curvature in 
an upward direction is revealed as the chiral limit is ap- 
proached, leading to the possibility of a larger infrared 
mass function for the lightest quark mass, despite the re- 
duction of the input bare quark mass. In contrast, a hint 
of downward curvature is observed for the most infrared 
points of the full QCD mass function as the chiral limit is 
approached. It is interesting that the nature of the curva- 
ture depends significantly on the chiral dynamics of the 
theory which are modified in making the quenched ap- 
proximation. Similar behavior is observed in the hadron 
mass spectrum where the coefficients of chiral nonana- 
lytic behavior can change sign in moving from quenched 
QCD to full QCD [Hil. 
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FIG. 3: Comparison of the unquenched (full QCD) and 
quenched quark propagator for non-zero quark mass. The 
mass function for the unquenched dynamical-fermion prop- 
agator has been interpolated so that it agrees with the 
quenched mass function for ma = 0.01 at the renormaliza- 
tion point, q — 3 GeV. For the unquenched propagator this 
corresponds to a bare quark mass of ma = 0.0087. 



V. CONCLUSIONS 

We have presented first results for the mass and wave- 
function renormalization functions of the quark propaga- 
tor in which the effects of 2 + 1 dynamical quark flavours 
are taken into account, fn contrast to the significant 
screening suppression of the gluon propagator in the in- 
frared [20|, the quark propagator is not strongly altered 
by sea quark effects. In particular, the renormalization 
function is insensitive to the light bare quark masses 
studied here, which range from 16 to 63 MeV, and also 
agrees well with previous quenched simulation results. 
Screening of dynamical mass generation in the infrared 
mass function is observed when comparing quenched and 
full QCD results at finite mass and in the chiral limit. 
The approach of the mass function to the chiral limit dis- 
plays interesting nontrivial curvature for low momenta, 
with the curvature in quenched and full QCD in opposite 
directions. 
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FIG. 4: Comparison of the unquenched (full QCD) and 
quenched quark propagator in the chiral limit. The renor- 
malization function is renormalized at q — 3 GeV. Whereas 
little difference is observed in the renormalization function, 
the mass functions indicates that dynamical mass generation 
is suppressed by the addition of quark loops. 
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APPENDIX: EXTRACTION OF THE SCALAR 
FUNCTIONS 

The Asqtad quark action is a staggered action us- 
ing three-link, five-link and seven-link staples as a kind 
of "fattening" to minimize quark flavor (often referred to 
as "taste") changing interactions. The three-link Naik 
term l'2Hl is included to improve rotational symmetry by 
improving the finite difference operator, and the five- link 
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FIG. 5: Renormalization (top) and mass (bottom) functions 
for four different quark masses in full QCD, where the va- 
lence and light-sea masses are matched. The wave-function 
renormalization function Z is renormalized at q = 3 GeV. 



Lepage term [24| is included to correct errors at low mo- 
menta that may be introduced by the above mentioned 
staples . The coefficients are tadpole improved and cho- 
sen to remove all tree-level 0(a 2 ) errors. 

At tree- level (i.e. no interactions, links set to the iden- 
tity) the staples in this action make no contribution, so 
the action reduces to the tree-level Naik action, 



S (0) = \ £x(*M*) \\{x{x + I*) - X(x /x)) 



X 

(A.l) 

where the staggered phases are: rj^(x) = (— 1)^ *" ' a and 



^(m) = [ 1 if v < M 

1 otherwise. 



(A.2) 
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FIG. 6: The chiral limit approach of the mass function for 
selected momenta. Results from quenched QCD simulations 
are illustrated at top whereas full dynamical-fermion QCD 
results are illustrated at bottom. Non-linear behavior is ob- 
served for the lowest momentum points, in opposite directions 
for quenched and full QCD. 



tion has the tree-level form 



? (or 

7 a/3 



(p; m) 

. , _ r9 1 

2j7 M )a/3 [g sin(p M ) - — sin(3p M ) 



mS. 



a/3 



= i £(7^)a/3sin(p Al )[l + isin 2 (p M )] +mS af3 . (A.3) 



where the a, /3 are themselves four- vectors: a M = 0,1, 
and likewise for f3; thus the quark propagator in Eq. IA.3I) 
is a 16 x 16 matrix. This familiar form is obtained by 
defining 



Safj — ILjJq^i mod 2 
(7 M ) Q/3 = (-1) Q, '^+C(m),/3- 



(A.4) 
(A.5) 



In momentum space, the quark propagator with this ac- The mod 2 in Eq. i|A.4|) ensures its validity in 



6 



Eq. ljA.5|) . The 7 satisfy 



f 7 1 * 

•y ' — ry - — *y — sy 
I fj, I (J, I fl I f-l ' 



(A.6) 
(A.7) 



forming a "staggered" Dirac algebra. 

Staggered actions are invariant under translations of 
2a, and the momentum on this blocked lattice is 

P^^f I m M = 0,...,^-l. (A.8) 

We calculate the quark propagator in coordinate space, 
G(x,y) = { X (x)x(y)}, (A.9) 

and obtain the quark propagator in momentum space 
by Fourier transform of G(x,0). To write the Fourier 
transform of the staggered field we write the momentum 
on the lattice 



27m, 



n M = 0,...,i p -l (A.10) 



so that k^ = p^ + ira^ and define J k = Then 



E> ( A - n ) 

X (x) = f e lk - x X (k) = / Y, ei{P+na) ' X X<*^ ( A ' 12 ) 



and 

G(x, y) = / exp {i(p + a7r) ■ x — i(l + [3ir) ■ j/} 

a/3 

x (Xa(p)X/j(0> (A.13) 
= E / e lp(:r ~ y) e l7r(Q - ;r ~' 9 - a) 5 Q;3 (p). (A.14) 

a/3 "'P 

Now it will be convenient to re-write this 

G{k) = G(l + ttS) = G s (l) = E e~ lkx G(x, 0) 

X 

= E / E exp { - ^ + 

a/3 J P x 

x exp{i(p + Tra)^}^^^) 
= E / d P i5asS a/ 3(p) (A. 15) 

(A.16) 



a/3 J P 

= E*wo- 



In the interacting case, the quark propagator asymp- 
totically approaches its tree-level value due to asymptotic 
freedom. At finite lattice spacing the actual behavior is 
closer to 



S(q;m) -> — S^(q;m/u ) 



(A.17) 

where uq is the tadpole (or mean-field) improvement fac- 
tor defined by 



u 



(<TrC/ plaq >fl\ 



(A.18) 



Assuming that the full lattice propagator retains its 
free form (in analogy to the continuum case) we write 

S J(P) = * E (%) a plM A (p) + B(p)6 a p (A.19) 

= Z- 1 (p) [i E (7 M ) Q/3 ^ &V) + M (P)^/?] . 

(A.20) 

where g is the tree-level momentum, Eq. ©. Combining 
this with Eq. I|A.15|I above, we can extract the scalar 
functions (which we now write in terms of g) as follows: 



G a (q)=Z(q)- 



-^(-^% + M(q) 



M 2 (q) 



(A.21) 



from which we obtain 



ETrG a (g) = 1W ( 



Z(q)M(q) 
q 2 + M 2 (q) 

lWJ3(q), 



(A.22) 



and 



zEE(- 1 )^^Tr[G Q (g)] = 167V c g 2 



Z(q) 



M 2 (q) 



= l6N c q 2 A(q). (A.23) 



Putting it all together we get 

A(q)=Z-\q) 
M(q) 



A(q) 



B(q) = 
M(q) = 



A 2 (q)q 2 +B 2 ( q ^ 
B(p) 



Z{q) A 2 (q)q 2 + B 2 (p) 

m 

A( q y 



(A.24) 
(A.25) 
(A.26) 



By calculating A, B instead of A, B, we avoid inverting 
the propagator. We calculate the ensemble average of A 
and B and thence M and Z. 
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